home *** CD-ROM | disk | FTP | other *** search
/ Tech Arsenal 1 / Tech Arsenal (Arsenal Computer).ISO / tek-04 / pgp23src.zip / SRC / MPILIB.H < prev    next >
C/C++ Source or Header  |  1993-05-09  |  17KB  |  469 lines
  1. /*    C include file for MPI multiprecision integer math routines.
  2.  
  3.     Boulder Software Engineering
  4.     3021 Eleventh Street
  5.     Boulder, CO 80304
  6.     (303) 541-0140
  7.  
  8.     (c) Copyright 1986-92 by Philip Zimmermann.  All rights reserved.
  9.     The author assumes no liability for damages resulting from the use
  10.     of this software, even if the damage results from defects in this
  11.     software.  No warranty is expressed or implied.
  12.  
  13.     These routines implement all of the multiprecision arithmetic necessary
  14.     for Rivest-Shamir-Adleman (RSA) public key cryptography, as well as
  15.     other number-theoretic algorithms such as ElGamal, Diffie-Hellman,
  16.     or Rabin.
  17.  
  18.     Although originally developed in Microsoft C for the IBM PC, this code
  19.     contains few machine dependencies.  It assumes 2's complement
  20.     arithmetic.  It can be adapted to 8-bit, 16-bit, or 32-bit machines,
  21.     lowbyte-highbyte order or highbyte-lowbyte order.  This version
  22.     has been converted to ANSI C.
  23.  
  24.     Modified 8 Apr 92 - HAJK - Implement new VAX/VMS primitive support.
  25.     Modified 29 Nov 92 - Thad Smith
  26. */
  27.  
  28. #include <string.h>
  29. #include "usuals.h"  /* typedefs for byte, word16, boolean, etc. */
  30. #include "platform.h" /* customization for different environments */
  31.  
  32. /* Platform customization:
  33.  * A version which runs on almost any computer can be implemented by
  34.  * defining PORTABLE and MPORTABLE, preferably as a command line
  35.  * parameter.  Faster versions can be generated by specifying specific
  36.  * parameters, such as size of unit and MULTUNIT, and by supplying some
  37.  * of the critical in assembly.  See the file platform.h for more 
  38.  * details on customization.
  39.  *
  40.  * The symbol HIGHFIRST, designating that integers and longs are stored
  41.  * with the most significant bit in the lowest address, should be defined
  42.  * on the command line for compiling all files, since it is used by files
  43.  * other than the mpilib routines.
  44.  */
  45.  
  46. #ifndef ALIGN
  47. #define ALIGN
  48. #endif
  49.  
  50. #ifndef PEASANT /* if not Russian peasant modulo multiply algorithm */
  51. #ifndef MERRITT /* if not Merritt's modmult */
  52. #ifndef UPTON   /* if not Upton's modmult */
  53. #ifndef SMITH
  54. #define SMITH   /* default: use Smith's modmult algorithm */
  55. #endif
  56. #endif
  57. #endif
  58. #endif
  59.  
  60. #ifdef SMITH
  61. #define UPTON_OR_SMITH  /* enable common code */
  62. #endif
  63. #ifdef UPTON
  64. #define UPTON_OR_SMITH  /* enable common code */
  65. #endif
  66.  
  67. #ifndef UNIT32
  68. #ifndef UNIT8
  69. #define UNIT16            /* default--use 16-bit units */
  70. #endif
  71. #endif
  72.  
  73. /***    CAUTION:  If your machine has an unusual word size that is not a
  74.     power of 2 (8, 16, 32, or 64) bits wide, then the macros here that
  75.     use the symbol "LOG_UNITSIZE" must be changed.
  76. ***/
  77.  
  78. #ifdef UNIT8
  79. typedef unsigned char unit;
  80. typedef signed char signedunit;
  81. #define UNITSIZE 8 /* number of bits in a unit */
  82. #define LOG_UNITSIZE 3
  83. #define uppermostbit ((unit) 0x80)
  84. #define BYTES_PER_UNIT 1 /* number of bytes in a unit */
  85. #define units2bits(n) ((n) << 3) /* fast multiply by UNITSIZE */
  86. #define units2bytes(n) (n)
  87. #define bits2units(n) (((n)+7) >> 3)
  88. #define bytes2units(n) (n)
  89. #endif
  90.  
  91. #ifdef UNIT16
  92. typedef word16 unit;
  93. typedef short signedunit;
  94. #define UNITSIZE 16 /* number of bits in a unit */
  95. #define LOG_UNITSIZE 4
  96. #define uppermostbit ((unit) 0x8000)
  97. #define BYTES_PER_UNIT 2 /* number of bytes in a unit */
  98. #define units2bits(n) ((n) << 4) /* fast multiply by UNITSIZE */
  99. #define units2bytes(n) ((n) << 1)
  100. #define bits2units(n) (((n)+15) >> 4)
  101. #define bytes2units(n) (((n)+1) >> 1)
  102. #endif
  103.  
  104. #ifdef UNIT32
  105. typedef word32 unit;
  106. typedef long signedunit;
  107. #define UNITSIZE 32 /* number of bits in a unit */
  108. #define LOG_UNITSIZE 5
  109. #define uppermostbit ((unit) 0x80000000L)
  110. #define BYTES_PER_UNIT 4 /* number of bytes in a unit */
  111. #define units2bits(n) ((n) << 5) /* fast multiply by UNITSIZE */
  112. #define units2bytes(n) ((n) << 2)
  113. #define bits2units(n) (((n)+31) >> 5)
  114. #define bytes2units(n) (((n)+3) >> 2)
  115. #endif
  116.  
  117. #define power_of_2(b) ((unit) 1 << (b)) /* computes power-of-2 bit masks */
  118. #define bits2bytes(n) (((n)+7) >> 3)
  119. /*    Some C compilers (like the ADSP2101) will not always collapse constant 
  120.     expressions at compile time if the expressions contain shift operators. */
  121. /* #define uppermostbit power_of_2(UNITSIZE-1) */
  122. /* #define UNITSIZE units2bits(1) */ /* number of bits in a unit */
  123. /* #define bytes2units(n) bits2units((n)<<3) */
  124. /* #define BYTES_PER_UNIT (UNITSIZE >> 3) */
  125. /* LOG_UNITSIZE is the log base 2 of UNITSIZE, ie: 4 for 16-bit units */
  126. /* #define units2bits(n) ((n) << LOG_UNITSIZE) */ /* fast multiply by UNITSIZE */
  127. /* #define units2bytes(n) ((n) << (LOG_UNITSIZE-3)) */
  128. /* #define bits2units(n) (((n)+(UNITSIZE-1)) >> LOG_UNITSIZE) */
  129. /* #define bytes2units(n) (((n)+(BYTES_PER_UNIT-1)) >> (LOG_UNITSIZE-3)) */
  130.  
  131. typedef unit *unitptr;
  132.  
  133.  
  134. /*--------------------- Byte ordering stuff -------------------*/
  135. #ifdef HIGHFIRST
  136.  
  137. /* these definitions assume MSB comes first */
  138. #define tohigher(n)        (-(n))  /* offset towards higher unit */
  139. #define pre_higherunit(r)    (--(r))
  140. #define pre_lowerunit(r)    (++(r))
  141. #define post_higherunit(r)    ((r)--)
  142. #define post_lowerunit(r)    ((r)++)
  143. #define bit_index(n)        (global_precision-bits2units((n)+1))
  144. #define lsbptr(r,prec)        ((r)+(prec)-1)
  145. #define make_lsbptr(r,prec)    (r) = lsbptr(r,prec)  
  146. #define unmake_lsbptr(r,prec) (r) = ((r)-(prec)+1) 
  147. #define msbptr(r,prec)        (r)
  148. #define make_msbptr(r,prec)    /* (r) = msbptr(r,prec) */
  149.  
  150. /* The macro rescale(r,current_precision,new_precision) rescales
  151.    a multiprecision integer by adjusting r and its precision to new values.
  152.    It can be used to reverse the effects of the normalize
  153.    routine given above.  See the comments in normalize concerning
  154.    Intel vs. Motorola LSB/MSB conventions.
  155.    WARNING:  You can only safely call rescale on registers that
  156.    you have previously normalized with the above normalize routine,
  157.    or are known to be big enough for the new precision.  You may
  158.    specify a new precision that is smaller than the current precision.
  159.    You must be careful not to specify a new_precision value that is
  160.    too big, or which adjusts the r pointer out of range.
  161. */
  162. #define rescale(r,currentp,newp) r -= ((newp) - (currentp))
  163.  
  164. /* The macro normalize(r,precision) "normalizes" a multiprecision integer
  165.    by adjusting r and precision to new values.  For Motorola-style processors
  166.    (MSB-first), r is a pointer to the MSB of the register, and must
  167.    be adjusted to point to the first nonzero unit.  For Intel/VAX-style
  168.    (LSB-first) processors, r is a  pointer to the LSB of the register,
  169.    and must be left unchanged.  The precision counter is always adjusted,
  170.    regardless of processor type.  In the case of precision = 0,
  171.    r becomes undefined.
  172. */
  173. #define normalize(r,prec) \
  174.   { prec = significance(r); r += (global_precision-(prec)); }
  175.  
  176. #else    /* LOWFIRST byte order */
  177.  
  178. /* these definitions assume LSB comes first */
  179. #define tohigher(n)        (n)     /* offset towards higher unit */
  180. #define pre_higherunit(r)    (++(r))
  181. #define pre_lowerunit(r)    (--(r))
  182. #define post_higherunit(r)    ((r)++)
  183. #define post_lowerunit(r)    ((r)--)
  184. #define bit_index(n)        (bits2units((n)+1)-1)
  185. #define lsbptr(r,prec)        (r)
  186. #define make_lsbptr(r,prec)    /* (r) = lsbptr(r,prec) */  
  187. #define unmake_lsbptr(r,prec) /* (r) = (r) */  
  188. #define msbptr(r,prec)        ((r)+(prec)-1)
  189. #define make_msbptr(r,prec)    (r) = msbptr(r,prec)
  190.  
  191. #define rescale(r,currentp,newp) /* nil statement */
  192. #define normalize(r,prec) prec = significance(r) 
  193.  
  194. #endif    /* LOWFIRST byte order */
  195. /*------------------ End byte ordering stuff -------------------*/
  196.  
  197. /*    Note that the address calculations require that lsbptr, msbptr, 
  198.     make_lsbptr, make_msbptr, mp_tstbit, mp_setbit, mp_clrbit, 
  199.     and bitptr all have unitptr arguments, not byte pointer arguments.    */
  200. #define bitptr(r,n) &((r)[bit_index(n)])
  201. #define bitmsk(n) power_of_2((n) & (UNITSIZE-1))
  202.     /* bitmsk() assumes UNITSIZE is a power of 2 */
  203. #define mp_tstbit(r,n) (*bitptr(r,n) &   bitmsk(n))
  204. #define mp_setbit(r,n) (*bitptr(r,n) |=  bitmsk(n))
  205. #define mp_clrbit(r,n) (*bitptr(r,n) &= ~bitmsk(n))
  206. #define msunit(r) (*msbptr(r,global_precision))
  207. #define lsunit(r) (*lsbptr(r,global_precision))
  208. /* #define mp_tstminus(r) ((msunit(r) & uppermostbit)!=0) */
  209. #define mp_tstminus(r) ((signedunit) msunit(r) < 0)
  210.  
  211.  
  212.     /* set working precision to specified number of bits. */
  213. #ifdef mp_setp
  214. void mp_setp(short nbits);
  215. #define set_precision(prec)    mp_setp(units2bits(global_precision=(prec)))
  216. #else
  217. #define set_precision(prec) (global_precision = (prec))
  218. #endif
  219.  
  220.  
  221. #ifdef PEASANT
  222.  
  223. /* Define C names for Russian peasant modmult primitives. */
  224. #define stage_modulus    stage_peasant_modulus
  225. #define mp_modmult    peasant_modmult
  226. #define modmult_burn    peasant_burn
  227. #define SLOP_BITS    PEASANT_SLOP_BITS
  228.  
  229. #else  /* not PEASANT */
  230. #ifdef MERRITT
  231. /* Define C names for Merritt's modmult primitives. */
  232. #define stage_modulus    stage_merritt_modulus
  233. #define mp_modmult    merritt_modmult
  234. #define modmult_burn    merritt_burn
  235. #define SLOP_BITS    MERRITT_SLOP_BITS
  236.  
  237. #else    /* not PEASANT, MERRITT */
  238. #ifdef UPTON
  239. /* Define C names for Upton's modmult primitives. */
  240. #define stage_modulus    stage_upton_modulus
  241. #define mp_modmult    upton_modmult
  242. #define modmult_burn    upton_burn
  243. #define SLOP_BITS    UPTON_SLOP_BITS
  244.  
  245. #else    /* not PEASANT, MERRITT, UPTON */
  246. #ifdef SMITH
  247. /* Define C names for Smith's modmult primitives. */
  248. #define stage_modulus    stage_smith_modulus
  249. #define mp_modmult    smith_modmult
  250. #define modmult_burn    smith_burn
  251. #define SLOP_BITS    SMITH_SLOP_BITS
  252.  
  253. #endif    /* SMITH */
  254. #endif    /* UPTON */
  255. #endif    /* MERRITT */
  256. #endif    /* PEASANT */
  257.  
  258.  
  259. #define mp_shift_left(r1) mp_rotate_left(r1,(boolean)0)
  260.     /* multiprecision shift left 1 bit */
  261.  
  262. #define mp_add(r1,r2) mp_addc(r1,r2,(boolean)0)
  263.     /* multiprecision add with no carry */
  264.  
  265. #define mp_sub(r1,r2) mp_subb(r1,r2,(boolean)0)
  266.     /* multiprecision subtract with no borrow */
  267.  
  268. #define mp_abs(r) (mp_tstminus(r) ? (mp_neg(r),TRUE) : FALSE)
  269.  
  270. #define msub(r,m) if (mp_compare(r,m) >= 0) mp_sub(r,m)
  271.     /* Prevents r from getting bigger than modulus m */
  272.  
  273. #define testeq(r,i)    \
  274.     ( (lsunit(r)==(i)) && (significance(r)<=1) )
  275.  
  276. #define testne(r,i)    \
  277.     ( (lsunit(r)!=(i)) || (significance(r)>1) )
  278.  
  279. #define testge(r,i)    \
  280.     ( (lsunit(r)>=(i)) || (significance(r)>1) )
  281.  
  282. #define testle(r,i)    \
  283.     ( (lsunit(r)<=(i)) && (significance(r)<=1) )
  284.  
  285. #define mp_square(r1,r2) mp_mult(r1,r2,r2)
  286.     /* Square r2, returning product in r1 */
  287.  
  288. #define mp_modsquare(r1,r2) mp_modmult(r1,r2,r2)
  289.     /* Square r2, returning modulo'ed product in r1 */
  290.  
  291. #define countbytes(r) ((countbits(r)+7)>>3)
  292.  
  293. /*    SLOP_BITS is how many "carry bits" to allow for intermediate
  294.     calculation results to exceed the size of the modulus.
  295.     It is used by modexp to give some overflow elbow room for
  296.     modmult to use to perform modulo operations with the modulus.
  297.     The number of slop bits required is determined by the modmult
  298.     algorithm.  The Russian peasant modmult algorithm only requires
  299.     1 slop bit, for example.  Note that if we use an external assembly
  300.     modmult routine, SLOP_BITS may be meaningless or may be defined in a
  301.     non-constant manner.
  302. */
  303. #define PEASANT_SLOP_BITS    1
  304. #define MERRITT_SLOP_BITS    UNITSIZE
  305. #define UPTON_SLOP_BITS    (UNITSIZE/2)
  306. #ifdef  mp_smul /* old version requires MS word = 0 */
  307. #define SMITH_SLOP_BITS    UNITSIZE
  308. #else           /* mp_smula or C version of mp_smul */
  309. #define SMITH_SLOP_BITS    0
  310. #endif /* mp_smul */
  311.  
  312. /*    MAX_BIT_PRECISION is upper limit that assembly primitives can handle.
  313.     It must be less than 32704 bits, or 4088 bytes.  It should be an
  314.     integer multiple of UNITSIZE*2.
  315. */
  316. #if (SLOP_BITS > 0)
  317. #define MAX_BIT_PRECISION (1280+(2*UNITSIZE))
  318. #else
  319. #define MAX_BIT_PRECISION 1280
  320. #endif
  321. #define MAX_BYTE_PRECISION (MAX_BIT_PRECISION/8)
  322. #define MAX_UNIT_PRECISION (MAX_BIT_PRECISION/UNITSIZE)
  323.  
  324.  
  325. /* global_precision is the unit precision last set by set_precision */
  326. extern short global_precision;
  327.  
  328.  
  329. /*    The "bit sniffer" macros all begin sniffing at the MSB.
  330.     They are used internally by all the various multiply, divide, 
  331.     modulo, exponentiation, and square root functions.
  332. */
  333. #define sniff_bit(bptr,bitmask)    (*(bptr) & bitmask)
  334.  
  335. #define init_bitsniffer(bptr,bitmask,prec,bits) \
  336. { normalize(bptr,prec); \
  337.   if (!prec) \
  338.     return(0); \
  339.   bits = units2bits(prec); \
  340.   make_msbptr(bptr,prec); bitmask = uppermostbit; \
  341.   while (!sniff_bit(bptr,bitmask)) \
  342.   { bitmask >>= 1; bits--; \
  343.   } \
  344. }
  345.  
  346. #define bump_bitsniffer(bptr,bitmask) \
  347. { if (!(bitmask >>= 1)) \
  348.   { bitmask = uppermostbit; \
  349.     post_lowerunit(bptr); \
  350.   } \
  351. }
  352.  
  353. /* bump_2bitsniffers is used internally by mp_udiv. */
  354. #define bump_2bitsniffers(bptr,bptr2,bitmask) \
  355. { if (!(bitmask >>= 1)) \
  356.   { bitmask = uppermostbit; \
  357.     post_lowerunit(bptr); \
  358.     post_lowerunit(bptr2); \
  359.   } \
  360. }
  361.  
  362. /* stuff_bit is used internally by mp_udiv and mp_sqrt. */
  363. #define stuff_bit(bptr,bitmask)    *(bptr) |= bitmask
  364.  
  365.  
  366. boolean mp_addc
  367.     (register unitptr r1,register unitptr r2,register boolean carry);
  368.     /* multiprecision add with carry r2 to r1, result in r1 */
  369.  
  370. boolean mp_subb
  371.     (register unitptr r1,register unitptr r2,register boolean borrow);
  372.     /* multiprecision subtract with borrow, r2 from r1, result in r1 */
  373.  
  374. boolean mp_rotate_left(register unitptr r1,register boolean carry);
  375.     /* multiprecision rotate left 1 bit with carry, result in r1. */
  376.  
  377. void mp_shift_right_bits(register unitptr r1,register short bits);
  378.     /* multiprecision shift right bits, result in r1. */
  379.  
  380. short mp_compare(register unitptr r1,register unitptr r2);
  381.     /* Compares registers *r1, *r2, and returns -1, 0, or 1 */
  382.  
  383. boolean mp_inc(register unitptr r);
  384.     /* Increment multiprecision integer r. */
  385.  
  386. boolean mp_dec(register unitptr r);
  387.     /* Decrement multiprecision integer r. */
  388.  
  389. void mp_neg(register unitptr r);
  390.     /* Compute 2's complement, the arithmetic negative, of r */
  391.  
  392. #ifndef mp_move
  393. #define mp_move(d,s)    memcpy((void*)(d), (void*)(s), \
  394.         units2bytes(global_precision))
  395. #endif  
  396. #ifndef unitfill0
  397. #define unitfill0(r,ct) memset((void*)(r), 0, units2bytes(ct))
  398. #endif
  399.  
  400. #ifndef mp_burn
  401. #define mp_burn(r) mp_init(r,0) /* for burning the evidence */
  402. #define mp_init0(r) mp_init(r,0)
  403. #endif
  404.  
  405. #define empty_array(r)  unitfill0(r, sizeof(r)/sizeof(r[0])/sizeof(unit))
  406.  
  407. void mp_init(register unitptr r, word16 value);
  408.     /* Init multiprecision register r with short value. */
  409.  
  410. short significance(register unitptr r);
  411.     /* Returns number of significant units in r */
  412.  
  413. int mp_udiv(register unitptr remainder,register unitptr quotient,
  414.     register unitptr dividend,register unitptr divisor);
  415.     /* Unsigned divide, treats both operands as positive. */
  416.  
  417. int mp_recip(register unitptr quotient,register unitptr divisor);
  418.     /* Compute reciprocal as 1/divisor.  Used by faster modmult. */
  419.  
  420. int mp_div(register unitptr remainder,register unitptr quotient,
  421.     register unitptr dividend,register unitptr divisor);
  422.     /* Signed divide, either or both operands may be negative. */
  423.  
  424. word16 mp_shortdiv(register unitptr quotient,
  425.     register unitptr dividend,register word16 divisor);
  426.     /* Returns short remainder of unsigned divide. */
  427.  
  428. int mp_mod(register unitptr remainder,
  429.     register unitptr dividend,register unitptr divisor);
  430.     /* Unsigned divide, treats both operands as positive. */
  431.  
  432. word16 mp_shortmod(register unitptr dividend,register word16 divisor);
  433.     /* Just returns short remainder of unsigned divide. */
  434.  
  435. int mp_mult(register unitptr prod,
  436.     register unitptr multiplicand,register unitptr multiplier);
  437.     /* Computes multiprecision prod = multiplicand * multiplier */
  438.  
  439. int countbits(unitptr r);
  440.     /* Returns number of significant bits in r. */
  441.  
  442. int stage_peasant_modulus(unitptr n);
  443. int stage_merritt_modulus(unitptr n);
  444. int stage_upton_modulus(unitptr n);
  445. int stage_smith_modulus(unitptr n);
  446.     /* Must pass modulus to stage_modulus before calling modmult. */
  447.  
  448. int peasant_modmult(register unitptr prod,
  449.     unitptr multiplicand,register unitptr multiplier);
  450. int merritt_modmult(register unitptr prod,
  451.     unitptr multiplicand,register unitptr multiplier);
  452. int upton_modmult(register unitptr prod,
  453.     unitptr multiplicand,register unitptr multiplier);
  454. int smith_modmult(register unitptr prod,
  455.     unitptr multiplicand,register unitptr multiplier);
  456.     /* Performs combined multiply/modulo operation, with global modulus */
  457.  
  458.  
  459.  
  460. int mp_modexp(register unitptr expout,register unitptr expin,
  461.     register unitptr exponent,register unitptr modulus);
  462.     /* Combined exponentiation/modulo algorithm. */
  463.  
  464. int mp_modexp_crt(unitptr expout, unitptr expin,
  465.     unitptr p, unitptr q, unitptr ep, unitptr eq, unitptr u);
  466.     /* exponentiation and modulo using Chinese Remainder Theorem */
  467.  
  468. /****************** end of MPI library ****************************/
  469.